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ABSTRACT 

Wc generalise the theory of mean-field galactic dynamos by allowing for temporal 
non- locality in the mean electromotive force (emf). This arises in random flows due 
to a finite response time of the mean emf to changes in the mean magnetic field, and 
leads to the telegraph equation for the mean field. The resulting dynamo model also 
includes the nonlinear dynamo effects arising from magnetic helicity balance. Within 
, <~| this framework, coherent large-scale magnetic spiral arms superimposed on the domi- 

Qh nant axially symmetric magnetic structure are considered. A non-axisymmetric forcing 

of the mean-field dynamo by a spiral pattern (either stationary or transient) is invoked, 
with the aim of explaining the phenomenon of magnetic arms. For a stationary dy- 
namo forcing by a rigidly rotating material spiral, we find corotating non-axisymmetric 
magnetic modes enslaved to the axisymmetric modes and strongly peaked around the 
corotation radius. For a forcing by transient material arms wound up by the galactic 
differential rotation, the magnetic spiral is able to adjust to the winding so that it 
resembles the material spiral at all times. There are profound effects associated with 
the temporal non-locality, i.e. finite 'dynamo relaxation time'. For the case of a rigidly 
rotating spiral, a finite relaxation time causes each magnetic arm to mostly lag the 
corresponding material arm with respect to the rotation. For a transient material spi- 
ral that winds up, the finite dynamo relaxation time leads to a large, negative (in the 
sense of the rotation) phase shift between the magnetic and material arms, similar to 
that observed in NGC 6946 and other galaxies. We confirm that sufficiently strong 
random seed fields can lead to global reversals of the regular field along the radius 
whose long-term survival depend on specific features of a given galaxy. 

Key words: magnetic fields - MHD - galaxies: magnetic fields - galaxies: spiral - 
galaxies: structure 



1 INTRODUCTION 

Nearby disc galaxies are known to typically have regular (or 
large-scale or mean) magnetic fields of 1-10 fiG in strength, 
which are coherent on the scale of the galaxies themselves 
(see Fletcher 2010; Beck 2012, for reviews). In many cases, 
galactic magnetic fields exhibit significant deviations from 
axial symmetry, in particular where the regular field is en- 
hanced in 'spiral magnetic arms' akin to the material arms. 
(By material arms, we specifically mean the regions where 
the densities of stars and gas are enhanced.) In most cases, 
there is clearly a relationship between the material spiral 
arms and regular magnetic field spiral arms. For instance, 
there exists a correspondence between the azimuthally av- 
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eraged pitch angle of the regular field and that of the mate- 
rial spiral (Fletcher 2010). This would be of little surprise, 
as the material arms naturally leave their imprint on the 
magnetic field, if not for the intriguing relation between 
magnetic and material arms, first discovered in the nearby 
galaxy NGC 6946 (Beck & Hoernes 1996). Here the mag- 
netic arms are located almost precisely between the mate- 
rial arms. They appear to be phase-shifted images of the 
material arms, with a negative phase shift in the sense of 
the galactic rotation (Frick et al. 2000; Beck 2007). Impor- 
tantly, the stronger tangling of the regular field by a (pre- 
sumably) more intense turbulence within the material arms 
cannot explain this phenomenon: both the total and regu- 
lar magnetic fields are stronger within the magnetic arms 
(Beck 2007). The nature of magnetic arms remains to be 
convincingly explained. 
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1.1 Magnetic arms and dynamo theory 

Turbulent mean-field dynamo theory has been successful 
in explaining the properties of the axisymmetric mode of 
regular fields in galaxies (Ruzmaikin et al. 1988, hereafter 
RSS88). 

Although growing non-axisymmetric modes can arise in 
this theory (beyond a certain radius in the disc where the 
rotational velocity shear is sufficiently small), they always 
have a lower growth rate than the axisymmetric mode if 
the galactic disc is axially symmetric (Baryshnikova et al. 
1987; Krasheninnikova et al. 1989). Thus, to explain non- 
axisymmetric modes as arising in an axisymmetric disc one 
must typically appeal to strongly non-axisymmetric seed 
fields and argue that the axisymmetric mode would not have 
had time to attain dominance. As we argue in this work, in- 
corporating field saturation through nonlinear feedback does 
not help the non-axisymmetric modes to become entrenched 
in an axisymmetric disc, and in fact causes them to eventu- 
ally decay. 

An alternative, of course, is to appeal to the devia- 
tions of the galactic discs from axial symmetry. Mestel 
& Subramanian (1991, hereafter MS91) and Subramanian 
& Mestel (1993, hereafter SM93) explored analytically and 
numerically the growth of non-axisymmetric modes under 
an enhancement of dynamo action along a spiral (presum- 
ably, but not necessarily, co-spatial with the material spiral) 
with a constant global pattern speed. Subsequently, Moss 
(1996, hereafter M96) (see also Moss 1998; Moss et al. 2001) 
carried out mean-field dynamo simulations to explore the 
effects of modulating various quantities, such as the a ef- 
fect, turbulent diffusivity, or the components of the mean 
velocity, along a spiral arm or bar. Some authors have also 
addressed directly the question of how the regular field could 
become enhanced in between the material arms (Shukurov 
2005; Rohde et al. 1999), though none of these explained 
the phenomenon of a substantial negative phase shift across 
a wide range of galactocentric distances, including those far 
away from the corotation radius, as reported by Frick et al. 
(2000). 

The earlier attempts to explain the phenomenon of 
magnetic arms encountered a fundamental difficulty of the 
same nature as in theory of spiral structure: it is feasible 
to maintain a significant phase difference between material 
and magnetic arms near the corotation radius, but the differ- 
ential rotation of the galaxy quickly winds up the magnetic 
arms to destroy the arms away from the corotation. In other 
words, in order to explain a constant phase shift between 
magnetic and material spiral arms, the former have to ro- 
tate rigidly. There are two important implications. Firstly, 
the large-scale magnetic field in magnetic arms cannot be 
frozen into the interstellar gas. Secondly, the governing equa- 
tions for the large-scale magnetic field must have wave-like 
solutions. The former implication is fully consistent with the 
mean-field dynamo theory where turbulent magnetic diffu- 
sivity can ensure the slippage of the large-scale magnetic 
field at appropriate time scales. However, the latter feature 
requires a significant extension of the existing theory since 
the leading growing modes of the standard mean-field dy- 
namo equations in a thin disc are non-oscillatory (in con- 
trast to spherical geometry, where dynamo waves can oc- 
cur). Thus, any attempt to explain magnetic arms within 



the framework of mean-field dynamos requires a substan- 
tial generalization of the theory. One such generalization is 
proposed here. 

1.2 Advances in dynamo and spiral structure 
theories 

Much of the recent work on mean-field dynamo theory has 
been focussed on the nonlinear regime, and the possible 
catastrophic quenching of the dynamo implied by magnetic 
helicity conservation. This has led to the development of the 
dynamical quenching theory (see the review of Branden- 
burg & Subramanian 2005a, hereafter BS05). In this the- 
ory, catastrophic quenching of the mean-field dynamo ac- 
tion is averted by a magnetic helicity flux, which transports 
small-scale magnetic helicity away from the region of dy- 
namo action. Dynamical quenching theory has been applied 
to local galactic dynamo models (Kleeorin et al. 2000; Vish- 
niac & Cho 2001; Kleeorin et al. 2002; Shukurov et al. 2006; 
Sur et al. 2007), and to axisymmetric discs (Kleeorin et al. 
2002; Smith 2012; Smith et al. 2012), but not yet to non- 
axisymmetric mean-field disc dynamos. Past work on the 
non-axisymmetric disc dynamos focused on the linear (kine- 
matic) regime, or relied on an approximate algebraic quench- 
ing formalism. 

Another recent development is the emergence of the 
minimal-r approximation (MTA) as a more general closure 
for mean-field electrodynamics that includes the quasilinear 
or first-order smoothing approximation (FOSA) as a limit- 
ing case (Kleeorin et al. 1996; Rogachevskii & Kleeorin 2000; 
Blackman & Field 2002; BS05). MTA is physically more 
appealing than FOSA because it takes into account the fi- 
nite response time of the mean electromotive force (emf) to 
changes in the mean magnetic field. This closure leads to 
new terms in the mean induction equation, which becomes 
a telegraph-type equation (Courant & Hilbert 1989), with 
second-order time derivative of the mean magnetic field. 
Separate considerations, motivated in part by the need to 
incorporate the non-locality in the dynamo coefficients, lead 
to essentially the same telegraph-type equation (Rheinhardt 
& Brandenburg 2012, see also Hughes & Proctor 2010). Such 
non-locality in time could lead to important effects in disc 
galaxies, where the product of the gas angular velocity io 
and correlation time of the turbulence r c may be of the or- 
der unity. Most notably, we find that it may lead to a large 
negative phase shift (in the sense of the galactic rotation) 
between the magnetic spiral arms and material spiral arms. 
For instance, for a spiral pattern rotating with pattern speed 
Q, one might expect a phase shift of order Qt c , due to the 
delay of order r c in the response of the mean emf to the 
enhanced a effect in the spiral arms. 

Yet another significant development has been the emer- 
gence of strong evidence that spiral patterns are not, at least 
in some cases, long-lived features rotating at a single con- 
stant pattern speed (Shetty et al. 2007; Dobbs et al. 2010; 
Sellwood 2011; Quillen et al. 2011; Roskar et al. 2011; Wada 
et al. 2011; Dobbs 2011; Kawata et al. 2011; Khoperskov 
et al. 2011). Theory, observations, and especially simula- 
tions of isolated and interacting galaxies point to a wide 
spectrum of possibilities, from relatively rigidly-rotating and 
long-lived patterns with fairly constant pattern speeds, to 
composite spirals comprised of multiple pattern speeds dom- 
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inating in different radial ranges, to material arms that ap- 
pear to rotate with the local gas velocity and thus quickly 
wind up. 

The goal of the present work is to draw together these 
recent developments in dynamo theory and spiral structure 
and apply the new ideas to examine non-axisymmetric regu- 
lar magnetic fields in disc galaxies. Here we focus on modes 
that are enslaved to the axisymmctric mode (and thus have 
the same growth rate in the kinematic regime); for a two- 
armed material spiral this would mean the quadrisymmet- 
ric mode (or m = 2 mode of the Fourier expansion) that 
corotates with the spiral pattern, as well as other even-m 
corotating modes, which are weaker. Here we discuss a nu- 
merical model, while Chamandy et al. (2012b, hereafter Pa- 
per II) presents a semi-analytical treatment of such modes. 
We leave the bisymmetric (m = 1) mode and other odd-m 
corotating modes to a forthcoming paper. Where the present 
work differs importantly from previous work is that here we: 

(i) incorporate MTA and explore the effects of a finite 
dynamo relaxation time; 

(ii) include dynamical quenching with an advective helic- 
ity flux for non-axisymmetric modes; and 

(iii) explore the effects of both rigidly rotating and tran- 
sient material arms on the mean magnetic field. 

In addition to this work on non-axisymmetric mean-field 
dynamos, we also briefly report on some new results from 
mean-field dynamo simulations which use an axisymmetric 
disc. 

The plan of the paper is as follows. We present in Sec- 
tion 2 the basic equations and introduce our mathematical 
approach. Wc outline the numerical model in Section 3. In 
Sect. 4 we discuss the findings from simulations of axisym- 
metric discs, while in Sections 5 and 6 we describe the re- 
sults of simulations of non-axisymmetric discs. A discussion 
of results and our conclusions are presented in Section 7. 



2 THE MEAN-FIELD DYNAMO AND THE 
TELEGRAPH EQUATION 

Consider the induction equation given by 

— = Vx (UxB-,,VxB), (1) 

where B is the magnetic field, U is the velocity field, r\ is 
the magnetic diffusivity. We follow the mean-field approach 
where the velocity and magnetic fields are each written as 
the sum of an average and a random component, 

B = B + b and U = U + u. (2) 

Here bar formally represents ensemble averaging but for 
practical purposes can be thought of as spatial average over 
scales larger than the turbulent scale but smaller than the 
system size. Substituting Eq. (2) into the induction equation 
leads to the mean-field induction equation (Moffatt 1978; 
Krause & Raedler 1980), 

r)R — — 

— = V x (U x B + £ - r,V x B) . (3) 

where 

£ = ulTb (4) 



is the mean electromotive force. 

Expressing £ in terms of the mean field is a standard 
closure problem. A relatively simple and widely used clo- 
sure is the quasilinear approximation, also known as the 
first-order smoothing approximation (FOSA; Moffatt 1978; 
Krause & Raedler 1980; BS05). In the quasilinear theory or 
FOSA, one neglects nonlinear terms in the evolution equa- 
tion for u and b, which, however, are retained if an evolution 
equation for £ is used instead (see below). For isotropic, he- 
lical turbulence, it leads to an expansion whose lowest-order 
terms are given by 

£ = aB - r,iV x B, (5) 

where, in the kinematic limit, a — o?k with 

a k = -|r c u ■ V x u, r) t = |r c u 2 , (6) 

and t c is the correlation time of the random flow. The tur- 
bulent transport coefficients a and r\ t are proportional, re- 
spectively, to the mean kinetic helicity and mean energy 
density of the turbulence. Below, we refer to the application 
of FOSA as the 'standard prescription'. 

An alternative treatment, suggested by Rogachevskii & 
Kleeorin (2000) and Blackman & Field (2002) was to replace 
the triple correlations which arise in the evolution equation 
for £ by a damping term proportional to £ itself (see BS05 
for a review). Under this approximation, called the minimal- 
r approximation (MTA) by BS05, one obtains, instead of 
(5), an evolution equation for £ , given by 

f = l(aB-*V*B)-f (7) 

where again in the kinematic limit, a and rj t are given by 
Eq. (6) and r is a relaxation time. For simplicity, this equa- 
tion has been derived assuming that r is scale-independent. 
When one takes into account the Lorentz force, the a- 
coefficient acquires an additional term proportional to the 
electric current helicity (Pouquet et al. 1976; Kleeorin & 
Ruzmaikin 1982; Gruzinov & Diamond 1994; Blackman & 
Field 2000; Radler et al. 2003; BS05), and then 

a = ak + a!m = -gTc[u-Vxu-b-Vxb]. (8) 

The T-approximation is motivated by the observa- 
tion that if, hypothetically, the mean fields were suddenly 
switched off then one would expect the mean emf to decay 
gradually, over a finite damping time r. This approximation 
has been tested in direct numerical simulations of forced tur- 
bulence (BS05; Brandenburg & Subramanian 2005b, 2007). 
These simulations of MTA find that r is positive and the 
associated Strouhal number is of order unity, ruofco ~ 1, 
where fco is the wavenumber corresponding to the correla- 
tion scale of the random flow, and uq is its rms velocity. In 
principle, the damping or relaxation time r can be different 
from the correlation time r c . For example, if r c is deter- 
mined by the frequency with which expanding supernova 
shocks encounter a given point in space then it could be 
shorter than r ~ (uofco) -1 (Shukurov 2004). Thus, we keep 
the ratio c T = r/r c as a dimensionless free parameter in the 
equations; for the numerical solutions we generally set it to 
unity. We then have 

(| + ±)* = ^(aB-*VxB). (9) 
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If the explicit time derivative is neglected (valid if £ varies 
on timescales long compared to the relaxation time r), and 
t is approximated as r c , then (9) reduces to the expression 
(5) obtained from the standard prescription. An alternative 
way of arriving at Eq. (9) (with c T — 1) is by keeping a 
time derivative of B in the expression (5) for £ in order to 
introduce non-locality in time (Rheinhardt & Brandenburg 
2012). 

We now apply the mean-field approach to the induction 
equation with the MTA closure. Operating on Eq. (3) with 
d/dt + 1/t, using Eq. (9), assuming U to be independent of 
time, and taking 77 and r\ t to be spatially uniform, we arrive 
at 

+ - [V x (U x B + c T aB) + (77 + cr-qt) V 2 B] . 

After multiplying through by r, this leaves us with an equa- 
tion containing new terms proportional to r (that do not 
emerge under the standard prescription) , 

d 2 B dB „ /== dB\ „ 2 dB 

+ V x (U x B + c T aB) + (77 + c T 77 t )V 2 B. 

(10) 

The same approach applied with the standard prescription 
gives the familiar result, 

^=Vx(UxB + aB) + ( r , + t,t)V 2 B, (11) 
at 

which is the r — > limit of Eq. (10) with c T — 1. Equa- 
tion (10) contains a second-order, as well as first-order time 
derivatives, and belongs to the class of equations known as 
the telegraph equation (e.g., Courant & Hilbcrt 1989). The 
second time derivative can lead to wave-like properties, with 
1/t as a damping coefficient. Unlike the mean induction 
equation (11), the telegraph equation (10) is not invariant 
under transformation to a rotating frame (see Paper II). 

3 NUMERICAL SOLUTIONS 

3.1 Method and approximations used 

We now examine the evolution of the mean magnetic field 
numerically. We solve Eq. (10) when r is finite, while 
Eq. (11) is solved for the case r — > 0. The components of 
Eq. (10) in cylindrical coordinates are given in Appendix Al. 
The saturation of the dynamo action is controlled by the he- 
licity conservation, outlined in Appendix A2, together with 
the governing equations. 

For convenience we define 

7= = V x£, 

and assume that 77, r; t , t and c T are constants. The quan- 
tity T has been introduced in addition to £ because this 
allows us to avoid having to impose boundary conditions 
on the components of £ and also to avoid applying the no- 
z approximation to the components of V x £ (see below). 
Alternatively, we could have solved for the vector potential 
A and £ , so the choice of using B and T is a matter of 
preference. 



We make the thin-disc approximation, which implies 
r~ 1 d/d(f> <§; d/dz and d/dr <C d/dz. We also assume that 
the variation of the rotational and radial velocities along 
the z-axis is negligible. Then we can neglect B z dU r /dz in 
the equation for B r and B z dU ^,/dz in the equation for B^. 
The thin-disc approximation and the condition V • B = 
together imply that B 2 /|B| -C 1. This approximation also 
allows us to ignore terms containing the (f>- and r-derivatives 
of aB z in Eqs. (A4) and (A5), and of B z in Eqs. (A8) and 
(A9), as well as the term £ Z B Z in Eq. (A7) for a m . Thus, the 
evolution equation for £ z is no longer needed, and all terms 
containing B z are eliminated from all the relevant evolution 
equations. Once equations for B r and B^ are solved, we can 
obtain B z from V • B = and T z from V • T — 0. As a 
bonus, we do not have to worry about satisfying V • B = 
and V • J- = numerically. 

Further, we write U$ = ruj, and take cj to be indepen- 
dent of 4>. We also adopt U r = 0. However, a vertical velocity 
U z can be essential for the dynamo action, so it is retained. 
The vertical flow advects magnetic field as well as magnetic 
helicity (Shukurov et al. 2006). Crucially, a magnetic helic- 
ity flux away from the dynamo region is required to alleviate 
the catastrophic quenching of the a effect (BS05). 

The no-2 approximation (Subramanian & Mestel 1993; 
Moss 1995; Phillips 2001) is used as discussed in detail in 
Appendix B, to approximate the ^-derivatives of the mean 
magnetic field. This reduces the three-dimensional problem 
to that in two dimensions (r and 0), thus greatly reducing 
the computational time required. The no-z approximation 
provides remarkably accurate solutions in one-dimensional 
dynamical quenching models (compare Shukurov et al. 2006 
and Sur et al. 2007) and there are no reasons to doubt that it 
performs equally well in higher-dimensional problems. Un- 
der this approximation, U z represents a vertically averaged, 
mass-weighted vertical velocity. 

3.2 Dimensionless governing equations 

We adopt the length scale and rms velocity of the largest tur- 
bulent eddies to be I — lo = 100 pc and u — uo — lOkrns -1 , 
respectively. These parameters are taken to be constant in 
space and time. This leads to the following estimate for the 
turbulent diffusivity: 

1 1 1 n26 2 —1 

f)t — 77t,o — 5 louo = 10 cm s . 

It is convenient to use dimensionless variables, with dis- 
tance along the z-axis measured in the unit of the charac- 
teristic half-thickness of the disc ho, horizontal lengths mea- 
sured in the radial disc scale length R, and time measured in 
typical vertical turbulent diffusion time £d,o = ^oAft.o- We 
take as fiducial values ho = 500 pc and R — 20kpc. For our 
purposes, it is not necessary to specify the magnetic field 
unit Bo, so we leave it arbitrary. Thus, we have lo = ho/ 5 
and Mo = 15/ioiao w i tn id ,o — 0.73 Gyr. 

Dimensionless variables, labelled with tilde, are thus in- 
troduced as h = hoh, ijt = r/tfiVt, r = Rr, t = ia,o£, t = 

td,0T, U = tJ^W, G — t^ Q G, Qk = fto^oSk, Ctm = hot^ am, 

U z = hot d ^ U z , B cq — BoB cq , B — BoB, £ = Bohot^t 
and T = Bot^T . The quantities ho and R enter only in 
the ratios A = ho/R and K — 2(ho/lo) 2 - The rotational 
shear is denoted G = rdui/dr, and B cq is a characteristic 
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field that may vary over space and time (normally associ- 
ated with the equipartition field; see below). Then, dropping 
tildes, our set of equations can be written in dimensionless 
form as 



dB r 

dt 

dt 
dT r 

dt 



dt 



dam 
dt 
d£ r 

~dt 

dt 



where 



dB r U z B r 
. G B r _ w ____ + ^ ) 

2C - R ^ R 

nh Ah 2 



+ c T r) t \ VB r y-^j- 



2C - R V R 

ixh Ah 2 



2 , ~ 2 9B, 
+ c T ri t \ 2 [PB 4> + — 



r 2 d(f> 



- J~6 



KrjtiS-B+TZjam) 



Ct m U z 



da n 



-r 1 {c T aB r - CrVtT^Bj, - £ r ^j 



7v I 3 —D 



V--— -— — — — 

r 2 r dr dr 2 r 2 d(f> 2 ' 



B r 



(12) 
(13) 

(14) 

(15) 
(16) 
(17) 

(18) 



3.3 Initial and boundary conditions 

The boundary conditions are B r = B^ = at r — and r = 
R, where the former condition applies as long as we consider 
non-axisymmetric modes that are localized far away from 
the rotation axis. We need not impose boundary conditions 
on the other variables as their spatial derivatives are not 
present in the equations. This includes the components of £ , 
since they are only evaluated for the purpose of calculating 
£ B, but not for calculating the components of Vxf . In the 
azimuthal coordinate, B r and B^, are subject to periodic 
boundary conditions, B r \<p=2TT = B r \<p=o and B^^-n = 

All variables vanish at t = except for B r and B^, 
whose seed values are chosen either to be Gaussian random 
fields for B r and zero for £>,*, or to have the functional form 



B 



r 
R 

B r = —B,),, 



r/R 



(So + Si COS (f>), 



(22) 
(23) 



where So and Si are dimensionless constants that control 
the amplitudes (relative to Bo) of the m = and m = 1 
components of the seed field. The solenoidality condition 
V ■ B = does not place any restriction on the form of the 
seed field, as B z can be obtained from V • B = 0. When we 
refer to the m = seed, this means So = 1 and Si = 0, 
while the m = 1 seed has So = and Si = 1. 

The code uses sixth-order finite differences for the spa- 
tial derivatives and a third-order Runge-Kutta routine for 
the time derivatives, using the same algorithms as the pub- 
licly available Pencil Code 1 . 



a = a k + a m , 



K = 2 



£ ■ B = S r B r + £<hB 



We have also neglected terms proportional to TZ^ 1 in 
Eqs. (12) and (13), as TZm 3> 1 in galaxies. We have retained 
such a term in Eq. (16); however, it can also be neglected if 
the flux is dominant, which is usually the case. 

For r — > 0, these equations reduce to the three standard 
equations of the slab dynamo: 

DBr 2C T -yj 7T 2 -yj 

Dt irh Ah 2 



+ c T ri t \ 



VB r 



2 dB 



DB^ - 2c T - n 2 - 

-DT= GB ^-^h aB ^- c ^^ B ' 

2 dB r 

VB 4> + 



Dam 
Dt 



+ C T 7J t \ 

■■-Kru 



r 2 d<f> 
c T a(B 2 r + B\) + IZ^an 
3V=D- 



(19) 



(20) 



(21) 



where 



_D 

Dt 



d d U z 
h to 1 

dt d<f> h 



3.4 The galaxy model 

For the galactic rotation curve, we use the Brandt curve, 



w(r) = 



wo 



[l + (r/r w ) a 



1/2 ■ 



(24) 



where u)o and r w are parameters, so that with this profile, 
u — > const as r — >■ (solid body rotation) and ua 1/r for 
r S> ruj (flat rotation curve). The rotational shear rate, 



duj 

G(r) = r- 



-u{r)- 



(r/r 



(25) 



' 1 + (r/r„) 2 ' 

tends to zero as r — > and G = —uj for r » r u , with the 
maximum magnitude of 2wo/(3v3) at r = y2r w . 

We model the disc half-thickness as a hyperboloid 
(RSS88), 

h(r) =hn [l + (r/r D ) 2 ] 1/2 , (26) 

where ho is the scale height at r = and ro controls the 
disc flaring rate. With this form, h — r ho = constant as 
r — > and h oc r for r >• to ■ 

We take the kinetic part of the a-coefficient to decrease 
with radius according to F. Krause's formula a k ~ l 2 co/h 
(RSS88; BS05). For the models of non-axisymmetric disc, 
we also impose a spiral profile on e*k. For a rigidly rotating 
spiral, we use 



Qk(r, <f>) — a(r) {1 + e a cos[n(4> — fit) — nr]} . 



1 http://pencil-code.googlecode.com (Brandenburg 2003) 



(27) 
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where the azimuthally averaged value of Qk, denoted with 
bar, is given by 

a(r) = l 2 u){r)/h{r), (28) 

e a sets the degree of deviation from axial symmetry, n is 
the number of spiral arms, fi is the angular velocity of the 
spiral pattern (i.e. the pattern speed), and k, negative for a 
trailing spiral, determines how tightly the arms are wound. 
For a spiral winding up with the differential rotation of the 
gas, we take 

a k (r, 4>) = a(r) [l + e a cos{n[<p + <j> - ut'}] , (29) 

where we have included an arbitrary phase </>o and defined 
t' = t — t on where t on is the time of onset of the a-spiral. We 
have also taken k — 0, so that the a-spiral actually starts 
off as a 'bar' at t' = 0. Wherever necessary, we also require 
that ak < u. 

The model also allows for a spiral modulation of the 
vertical advective velocity U z of a similar form: 

Uz = Uz,o{l + e v cos[n(<j) - Qt) - nr]}, (30) 

were U Zl o is the azimuthally averaged value of U z . 

Following other authors (e.g. RSS88) we define the tur- 
bulent magnetic Reynolds numbers to quantify the strengths 
of the differential rotation, a effect and vertical velocity: 

Ru = Gh 2 /r] t , R a =a k h/rj t , Ru = U z h/rj t . 

The local dynamo number, 

D = R a Ru = a k Gh s /ri 2 , 

is a dimensionless measure of the intensity of the (aw) dy- 
namo action at a given radius. For the mean field to grow 
due to the local acj-dynamo action alone, it is required that 
D < D cr « -10 (RSS88, see also Paper II). In the limit 
r >• r w , an axisymmetric disc has D ~ —9(luIi/u) 2 . 

The functional form of the characteristic magnetic field 
is adopted as 

B ca = B e- r/R ; (31) 

if appropriate, this can be identified with the magnetic field 
strength corresponding to the energy equipartition with the 
turbulence. Nonlinear dynamo effects are expected to be- 
come pronounced as soon as |B| ~ B cq . 

3.5 Models explored and the representation of 
the results 

The various numerical models considered here are summa- 
rized in Table 1. Many parameters are common to all models. 
For the disc scale length, we choose ru = R/2 = lOkpc and 
set its scale height to be h — ho = 500 pc at r = R/2 = 
lOkpc, so that /id = ho/\[2 = 354 pc. 

We adopt r w = 0.1R = 2kpc and set the circular speed 
to U4, = 375MJ,o = 250 km s -1 at r = R/2 = lOkpc, 
which gives ojo = 96£7q = 130 kms -1 kpc -1 . This implies 
that uj = 18.75^o = 25kms~ 1 kpc~\ G ~ — 18*^" = 
-24kms~ 1 kpc -1 , a = 0.75hot^ = 0.5 kms -1 and the 
azimuthally averaged dynamo number D ~ —13.5 at r = 
R/2 = lOkpc. 

The resulting radial profiles of h, _R„, R a , and D are 
shown in Fig. 1. 



Table 1. List of numerical models. Models A-D have an axisym- 
metric disc. The value of Ru given is its value at r = R/2 = 
lOkpc, and, for Model F, its azimuthal average. 'R' in Column 2 
indicates that the seed is Gaussian random noise in B r , whereas 
numbers and 1 refer to the azimuthal symmetry of the seed 
field, m = and m = 1, respectively. In Column 7, 'T' indicates 
that the imposed spiral is transient, while 'W indicates that the 
a-spiral begins as a bar and winds up with the gas. In all models, 
runs were performed (i) without r (r — ► limit), and (ii) with 
t = Iq/uo. All models used R = 40ho = 20kpc, as well as the 
profiles for the disc thickness and gas rotation velocity discussed 
in Section 3.4. 
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0.45 
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K 
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0.45 


0.5 
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L 





0.2 


0.45 


0.5 





20 


4 


M 





0.2 


0.45 


0.5 





20,T 


2 


N 





0.2 


0.45 


0.5 





T,W 


2 



The corotation radius for models with steady spiral forc- 
ing is chosen as r COI = 8kpc, which is close to the observa- 
tional estimates for the Milky Way (Gerhard 2011; Acharova 
et al. 2011). This leads to a pattern speed fl = 23£Jq = 
31 kms -1 kpc -1 . Moreover, a reasonable estimate for r is 
given by r = lo/uo = (Zo/3/io)^d,o = 0.013£a,o- It is possible 
that r is much smaller than our simple estimate for some 
galaxies and much larger for others. Therefore, we also con- 
sider the t — y case, as well as r = 2lo/uo and r = 4/o/wo 
for some models. 

We decompose each component of the mean magnetic 
field B into a cosine Fourier series, with certain phases (j)^ 
(see below), 

00 

Bi(r,4>,t) = J2 Bi m) (r,t)cos{m[(^-^(r,t)]}, (32) 
where i — r,4>, Thus, 

B i 0) = ^j o Bi{r,<j>,t)d4>, (33) 

and, for m > 0, 

B { r ] = Bi(r, 4>, t) cos {m [0 - ^ (r, *)] } #. 

(34) 

The phase (j)^ is obtained by trying all possible values 

and choosing the one which maximizes B\ . Using their 
phases, we also determine the rotation rates of the various 
Fourier modes by numerically differentiating the phase with 
respect to time. 

The azimuthal average of the magnetic energy can be 
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Figure 1. Inputs used for the basic non-axisymmetric Model E. 
(a) Disc scale height h, given by Eq. (26). (b) The dimensionless 
quantity — R^ = Gh 2 /r/ t with G given by Eq. (25). (c) The di- 
mensionless quantity R a = a-^h/rjt, with given by Eq. (27), at 
(j> = (solid), azimuthal maximum/minimum (dotted), and az- 
imuthal mean (dash-dotted), (d) Negative of the dynamo number 
—D = —RaRu at (j> = (solid), along with its azimuthal extrema 
(dotted) and mean (dash-dotted). The azimuthal mean is used in 
axisymmctric models A-C. Negative of the approximate critical 
dynamo number —D cr ~ 10 is shown dashed for reference. 



written as 



where 



E 



(0), 



E(r,t)= ^£ (m) M), 
,t) = l{[Bi°Hr,t)} 2 + [B$\r,t)]- 



(35) 



(36) 



and, for m > 0, 



E 



(m) 



(r, t) = ^ { [B^ (r, «)] 2 + [ B <?> (r, t)] ' } , (37) 

where the additional factor 1 /2 in the latter equation arises 
from averaging cos 2 [m((j) — (j)^)] over the interval (0, 2tt). 
Averaging over the area of the disc, we obtain the average 
magnetic energy in mode m, 

r-R 



(E (m) )(t) = 



B? 



E (m) (r,t)rdr. 



(38) 



The numerical mesh used has n r — 200 nodes in radius 
and n$ = 120 in azimuth for most runs, though for some 
runs we used n r — 100 and = 60. 



4 DYNAMO IN AN AXISYMMETRIC DISC 

A number of interesting questions can be addressed by con- 
sidering dynamo action in an axisymmetric disc. One of 
them is how a whole young galaxy becomes magnetized with 



a coherent field given that the radial diffusion time across a 
galaxy, tn ~ R 2 /r/ t — 3x 10 2 Gyr, is typically larger than the 
age of the Universe. We examine this in the context of the 
dynamical quenching model of the galactic dynamo, using 
the no-z approximation. 



4.1 Axisymmetric solutions 

4-1.1 The role of the helicity flux 

In local dynamo models (one-dimensional, retaining the z- 
dependence of the solution alone or zero-dimensional, us- 
ing the no-z approximation) advection of magnetic helicity 
through the disc surface has been shown to alleviate the 
catastrophic nonlinear quenching of the dynamo (Shukurov 
et al. 2006; Sur et al. 2007). In the absence of any helic- 
ity flux, the mean magnetic field does grow to a fraction 
of the equipartition value, until a m cancels a^, after which 
magnetic field decays. 

The result which obtains when the radial dimension is 
included is more interesting, even in the absence of a helicity 
flux. This case is studied using Model A (using the no-z ap- 
proximation), with the evolution of the radial profile of the 
mean field strength shown in Fig. 2 for the case r — > 0. No 
notable differences exist between the r — > and r = Iq/uq 
cases. For both cases, the kinematic global growth rate is 
virtually the same, T = 5.7£7q ~ 7.8 Gyr -1 . Since the local 
growth rate of the field depends on the dynamo number, 
which in turn depends on r, the radius, the field maximum 
travels in radius. More specifically, the maximum is first lo- 
calised where the dynamo number is the largest (jm, say), 
but, at later time, at radii tm — Ari and tm + Ar2, etc. 
This leads to two rings of enhanced mean field, one mov- 
ing inward and the other moving outward (the outer ring is 
more prominent in Model A). Because of the catastrophic 
a-quenching, the mean field eventually becomes negligible. 
(See Moss et al. (1998) for a discussion of propagating mag- 
netic fronts in disc galaxies.) 

On the other hand, if the vertical advection, and thus 
helicity flux, is above a certain threshold (but not too large), 
catastrophic quenching is averted. Then magnetic field can 
persist at about the same strength in a wide radial range, 
gradually spreading out until it occupies the entire region 
of the disc where the local dynamo number is supercritical. 
This is the situation illustrated with Model B, as shown in 
Fig. 3 (r -> case), for R v = 0.45 at r = R/2 = 10 kpc 
and U z — 0.3 km s" 1 at all r. In this model, both r — > and 
r = lo/uo cases are again very similar, and the kinematic 
global growth rate for both cases is F = 5.1tJ ~ 7.0 Gyr -1 . 
This is slightly smaller than for Model A, since the vertical 
flux removes large-scale magnetic field, rendering the dy- 
namo less efficient in the kinematic stage. This model clar- 
ifies how the whole disc becomes magnetized. It does so by 
first reaching significant strengths in regions where the dy- 
namo number is largest. Growth then saturates at each ra- 
dius at a fraction of the local equipartition value, provided 
the vertical velocity (or helicity flux) is favourable for the 
dynamo action. The role of the radial diffusion is merely 
to couple the dynamo action at different radii as to lead to 
a magnetic structure growing at a single rate. Importantly, 
the entire disc becomes magnetized over a timescale much 
shorter than the radial diffusion time. 



© 2012 RAS, MNRAS 000, 1-25 



8 L. Chamandy, K. Subramanian and A. Shukurov 



0.5 F 



0. 1 



0.3 - 



=3 



0.2 - 



0.1 -: 



0.0E 



t = 3.1' Gyr 

t = 3.5 Gyr 

* = 4.2 Gyr 

t = 4.6 Gyr 

t = 5.5 Gyr 

t = 7.0 Gyr 

J = 10.0 Gyr 
0-.5fl e /D a 




12 



«5 

"3- 



0.5 



0.4 r 



0.3 



0.2 



0.1 



0.0 



-0.1 



-0.2 



I \ 
I \ 



•— * = (B r /B„ is plotted): 

- t = 0.7 Gyr 
■ ■- t = 5.8 Gyr 

— t = 11.0 Gyr 
' ~ = 11.0 Gyr (m=0 seed): 



2 



4 6 8 10 12 
r (kpc) 



Figure 2. Evolution of the field strength B = (B% + B^) 1/2 in 
Model A (axisymmetric disc, with Rjj = 0). The mean field is 
confined to radially expanding and contracting rings. 
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Figure 3. Evolution of the field strength in Model B (axisymmet- 
ric disc, with Ru = 0.45 = 0.3kms _1 ). The mean field spreads 
radially to fill the disc. 

4.2 Reversals of the magnetic field 

If the seed magnetic field is weak enough, the steady- 
state magnetic configuration is independent of its form and 
strength. However, a relatively strong initial magnetic field 
can affect the steady-state magnetic configuration if non- 
linear dynamo effects become important before the leading 
dynamo eigenfunction (normally represented by an axisym- 
metric magnetic field without any reversals along the radius) 



Figure 4. Evolution of the magnetic field that starts as a random 
seed magnetic field (Model C) in the r — > case. Curves show 
at the azimuth (f> = at various times, with the exception of the 
field at t = 0, where B r is shown (since B^ = 0). Otherwise, B r is 
omitted for the sake of clarity, but vanishes at the same set of radii 
as B^, indicating reversals of the field at these locations. These 
reversals can survive for several Gyr, slowly moving outwards. For 
the final time, the solution obtained with an m = seed is shown 
as a thin dotted line for reference. 

can become dominant (Shukurov 2005; Moss et al. 2012). In 
particular, a random seed magnetic field can lead to long- 
lived reversals in the (quasi-)steady state magnetic configu- 
ration, either global (Poezd et al. 1993) or localised in both 
radius and azimuth (Bykov et al. 1997). A suitable random 
seed magnetic field can be readily provided by the fluctu- 
ation dynamo action (Poezd et al. 1993). Earlier results in 
this area have been obtained with a heuristic algebraic non- 
linearity in the mean-field dynamo equation. Here we revisit 
this idea, but now with the physically motivated dynamic 
nonlinearity and finite r. 

In Model C, we have chosen a random seed specified as a 
two-dimensional Gaussian random field with the root-mean- 
square value of about one tenth the equipartition value (i.e., 
a fraction of fj,G in the Solar neighbourhood) 2 We also tried 
ten times larger seeds and found almost identical results. 
The number of radial grid points in this models is n r = 100, 
corresponding to a resolution 0.2 kpc in radius, comparable 
to the correlation scale of the interstellar turbulence. The 
results are shown in Fig. 4. Firstly, large scale fields can de- 
velop over kpc scale regions, even in the outer disc, on Gyr 
timescales. More interestingly, it can be seen that reversals 
develop and persist for several Gyr in the nonlinear regime. 
(Because of the random nature of the problem, other ran- 
dom seeds of the same strength or different resolutions can 
result in field configurations without reversals.) The rever- 

2 The solenoidality of the seed field is ensured by an appropriate 
choice of B z . 
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Figure 5. The strength of the mean magnetic field (colour coded) at (a) t = 0.7 Gyr, (b) 2.8 Gyr, (c) 7.3 Gyr, and (d) 11.0 Gyr in 
Model C (random seed field) for the case r — > 0. The colour of the central region has been saturated to enhance visibility of the magnetic 
structure in the outer disc. Field vectors are also shown, with tail length proportional to the magnitude of magnetic field, for r > 5 kpc 
only, to avoid clutter. Global reversals whose positions evolve with time are visible (for example, the central dark ring which moves out 
with time). 



sals are also apparent in a time sequence of the magnetic 
field of Model C shown in Fig. 5. At least one reversal in 
the regular magnetic field has been observed in the Milky 
Way (Van Eck et al. 2011). It can be seen from both figures 
that as time goes on, the reversals become global in nature, 
propagate outward on a timescale of several Gyr, and in the 
process decrease in number. By t = 6£d,o = 4.4 Gyr all rever- 
sals inside r = 16 kpc are global in nature (i.e. occur at all 
azimuth for a given radius). In the r = lo/uo case, the noise 
from the random seed takes much longer to dissipate, and 
global reversals are only apparent much later on, at large 
radii where the field is very weak. The occurrence of the 
reversals and their sensitivity to model parameters need to 
be explored further. In particular, Poezd et al. (1993) found 
that the persistence of the reversals at the galactic lifetime 
scale strongly depends on the rotation curve. Our experi- 
ments reported here suggest that memory effects (finite r) 
are also important. 

4.3 Non-axisymmetric magnetic fields in an 
axisymmetric disc 

To test the performance of our code with non-axisymmetric 
solutions, we have repeated several of the 2D kinematic- 
regime simulations of M96 with the same resolution (41 grid- 
points in r and 61 in <j>), the same parameters, R = 20/io, 
cjot^ 1 = 10, r w = 0.2R, R a — const = 1, and without either 
disc flaring (h = const) or radial variation of the a effect, 
or the extra factors oi 2/n and 7r 2 /4 used to refine the no-z 
approximation elsewhere in this paper. M96 only considered 
axisymmetric and bisymmetric modes, m = 0, 1, to obtain 
their growth rates To = O.SOt^" 1 and Fi = O.lbt^ 1 , respec- 
tively; the bisymmetric mode is localised near the edge of 
the disc and rotates at the rate Q\ — 2A6t~[ . We con- 
firm the existence of these modes and find good agreement 
with the results of M96. Our results are comfortably simi- 
lar: To = 0.8WJ 1 , Ti = 0.15^ and fii = 2.53i~ J ; the slight 
differences are likely due to the different algorithms used. 

We note, however, that some conclusions of M96 may 
rely on the details of an oversimplified model that was used. 
In particular, the m = 1 mode grows because is in- 



dependent of r rather than being proportional to uj as in 
Krause's law. This allows the a effect to dominate over the 
rotational shear at large r, which favours the growth of non- 
axisymmetric modes as in a 2 versus ato dynamos (Raedler 
1986). The factors 2/ir and 7r 2 /4 in the no-z approximation, 
absent in that model, reduce the a effect and enhance the 
importance of the turbulent diffusion, which also has an ad- 
verse effect on the growth of non-axisymmetric modes. We 
find that the non-axisymmetric modes decay in a more re- 
alistic model used here and for the parameters used in most 
of our numerical models (see Sect. 3.5). We also find, how- 
ever, that the m = 1 mode can grow if we adopt somewhat 
modified parameters, e.g., a smaller rms turbulent velocity 
u — 5kms _1 and a larger turbulent scale length I = 150 pc. 
This makes sense because lowering u weakens the radial dif- 
fusion while increasing I enhances R a , without affecting the 
shear. With this new set of parameters, has to be trun- 
cated in the central region to ensure that a ^ |«. The effect 
of the truncation is to shift the m = 1 mode to larger radius 
(near the peak of R a ) than it otherwise would be. It is worth 
noting that higher m modes are even more difficult to excite 
in an axisymmetric disc and would require the parameters 
to be more extreme still. 



4-3.1 Nonlinear competition of non-axisymmetric modes 

In Model D, we focus on the case r — > 0. The m = 1 mode is 
readily excited with the kinematic growth rate Ti = 1.7t^ . 
The growth rate for the axisymmetric mode for these param- 
eters is Fo = 4.7t^ q. Since Fq > Ti with a significant margin, 
magnetic field in a mature dynamo hosted by an axisymmet- 
ric disc will be axisymmetric unless the seed magnetic field 
is strongly non-axisymmetric. The latter is, in fact, quite 
plausible if the seed field arises, for example, from a puta- 
tive intergalactic field captured as the disc galaxy forms. In 
order to address this point, we have run Model D, starting 
from a purely m = 1 seed, with the m = mode being 
seeded only by the numerical noise, as in M96. We adjust 
the strength of the seed such that nonlinear dynamo effects 
become significant at various important stages: (i) when the 
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Figure 6. Evolution of the normalised magnetic field strength 
averaged over the disc, in the m = (thick) and m = 1 (thin) 
modes, for a purely bisymmetric (m = 1) seed magnetic field 
(Model D). Different line styles show three cases with identical 
parameters, except for the magnitude of the seed field, so that 
the dynamo action saturates (i) when the m = mode dominates 
(solid), (ii) when m = 1 dominates (dashed), and (iii) right from 
t = (dash-dotted). Only in case (ii) does a significant m = 1 
component survive for several Gyr. 

m = mode has already come to dominate, (ii) when it is 
still weaker than the m = 1 mode, and (iii) right from t — 0. 

The results are shown in Fig. 6. Thick lines correspond 
to the normalised magnetic field strength (averaged over the 
area of the disc) in the m = mode while thin lines corre- 
spond to that in the m = 1 mode. In Case (i), shown by the 
solid lines, the saturation of the (essentially axisymmetric) 
field clearly causes the m = 1 mode to decay since its growth 
requires a stronger (unquenched) dynamo action. 

In Case (ii) (dashed lines), the early saturation of the 
(stronger) m = 1 mode still does not prevent the m — 
mode from growing and evolving as in Case (i). This hap- 
pens because the magnetic field in between the extrema of 
the non-axisymmetric magnetic field remains relatively weak 
when the m = 1 mode ceases to grow. Therefore, the grow- 
ing m — mode is supported by the local dynamo action 
in those regions, to fill the gaps until the field is saturated 
everywhere to become nearly axisymmetric. However, the 
m = 1 mode, though subdominant, does not immediately 
decay but remains strong for several Gyr. This feature can 
be understood from Fig. 7, which shows snapshots of B<t, / Bo 
at various times. Initially, the m — 1 mode is growing expo- 
nentially. Like the axisymmetric mode in Case (i), the m = 1 
mode begins to spread radially as soon as nonlinear effects 
become important. At any given radius, the m — mode 
eventually catches up and begins to dominate there. Thus, 
the outward spreading of the dominant m — mode is pre- 
ceded by the outward spreading of the m = 1 mode. As a 
result, a ring of m = 1 is formed (blue and red spirals in 
Fig. 7), surrounding a central disc of m = 0, both expand- 
ing outward with time. Eventually, the axisymmetric mode 
comes to dominate everywhere and the m = 1 decays more 
quickly (for t > 12 Gyr). 

In Case (iii) (dash-dotted lines), we find that the m = 1 
mode first decays and then remains much weaker than the 
m = mode. 

As far as the overall survival of non-axisymmetric 
modes in an axisymmetric disc is concerned, we can con- 
clude that: 



(i) a rather special parameter combination is needed for 
a non-axisymmetric mode to grow; 

(ii) for a non-axisymmetric mode to persist for several 
Gyr in the nonlinear stage, it must reach the saturation 
strength before the m = can do so; given that non- 
axisymmetric modes have smaller growth rates, this implies 
that they must be much stronger initially; 

(iii) in any case, all non-axisymmetric modes eventually 
decay in an axisymmetric disc. 

These conclusions make it all the more reasonable to sug- 
gest that deviations from perfect axial symmetry of the un- 
derlying disc are required to explain the prevalence of non- 
axisymmetric regular fields in many galaxies. It is the mech- 
anism of 'spiral forcing' of the dynamo to which we now 
turn. 



5 MAGNETIC ARMS ENSLAVED TO A 
STATIONARY SPIRAL PATTERN 

In this section, we consider the evolution of the mean mag- 
netic field in a non-axisymmetric disc where the a effect 
is modulated by a stationary spiral pattern. As above, our 
analysis includes the dynamic nonlinearity, thus extending 
fully into the saturated states of the dynamo action, and al- 
lows for a finite dynamo relaxation time r. The results at the 
kinematic stage can be compared with the semi-analytical 
solutions of Paper II. 

As a consistency check, we have first made sure that 
the standard dynamo solutions, obtained by solving Eq. (11) 
are obtained when solving Eq. (10) with r approaching zero. 
Indeed, we found good agreement with them for r = 1CF 3 . 
To ensure that the results are not sensitive to the location of 
the outer boundary, we compared solutions with outer radius 
of R = 15kpc and 30kpc, with the standard value R = 
20kpc (with the radial resolution modified proportionately). 
The results are not sensitive to this adjustment. 

Our standard model is labelled F, with parameters 
given in Table 1. To explore the parameter space, we have 
also considered Models F-L. 

5.1 Kinematic solutions 

The evolution of the normalized magnetic field strength in 
various even-m modes, averaged over the area of the disc, 
is shown, for Model E, in Fig. 8. In the kinematic regime, 
all modes have about the same growth rate T ~ 5.1tJ ~ 
7.0 Gyr -1 for both r — > and r = Iq/uq cases. This growth 
rate is also in close agreement with the growth rate for the 
axisymmetric case (Model B). This is not surprising be- 
cause the average field in the disc is used to calculate the 
growth rate, and the field is dominated by the contribution 
at r <C r CO r, where non-axisymmetric modes are negligible. 
The even-m modes grow together with the m — mode be- 
cause they are driven by (enslaved to) it. The m = 2 and 
m = 4 modes corotate with the spiral pattern. The odd 
modes are so weak that they are dominated by numerical 
noise, and thus also grow with the m — mode. 

We illustrate in Fig. 9 various aspects of the solution at 
3td,o after the simulation has begun (or t = —5.8 Gyr in the 
plots), when the dynamo is still in the kinematic regime, and 
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Figure 7. The </>-component of the magnetic field in Model D, i.e. an m = 1 seed field in an axisymmetric disc, at (a) t = 2.2 Gyr, (b) 
3.7 Gyr, (c) 5.1 Gyr, and (d) 6.6 Gyr. The results shown are for Case (ii), as described in the text, with the initial average strength of 
the m = 1 mode of about 10 — 4 i?o (dashed curves in Fig. 6). The bisymmetric (m = 1) mode survives for several Gyr in the nonlinear 
regime in the form of an expanding annulus surrounding an expanding m = disc. 



axisymmetric part, 




Figure 8. Magnetic field strength in the m = (solid), m = 2 
(dashed) and m = 4 (dash-dotted) modes of the mean magnetic 
field, averaged over the entire disc and normalized to Bo, the 
equipartition field strength at r = 0, for Model E. Results ob- 
tained for t = lo/uo are shown in thick black, while those for 
t — y are in thin red. All modes have exponential kinematic 
growth rate r< m ) = 7.0 Gyr" 1 . For the convenience of presenta- 
tion, time has been rescaled so that the simulation in fact starts 
from t = —8 Gyr (not shown). 



one can compare these results with the semi-analytical solu- 
tions of Paper II. In Fig. 9a, B^, (solid) and —B r (dashed) 
have been plotted for both r — > (red) and r = Iq/uo 
(black) cases, at 4> = <^w, one of the two azimuthal angles 
where the a-spiral crosses the corotation circle. As expected, 
the field is strongest inside the corotation radius, where the 
magnitude of the dynamo number is largest - this is clear 
from Fig. 9a. Importantly, however, there is an excess of 
magnetic energy near the corotation radius, as compared to 
what is obtained in an axisymmetric disc, mainly due to the 
presence of strong additional m = and m — 2 components 
near r = r cor . The m = 4 component of the field is found to 
always be much weaker than m = 2. 

The presence of the non-axisymmetric component near 
corotation can be more clearly seen in Fig. 9b. There we 
plot the ratio of the non-axisymmetric part of B^ to the 



8 = 



Bj, — B 



B 



(()) 



at the azimuth 4> — <f> COI . For both vanishing and finite r, 
the eigenfunction is important within an annulus of width 
m 4 kpc, centred near the corotation circle. Note that in 
both cases, the radial phase varies much more rapidly than 
that of a (plotted as a thin dash-dotted line in the figure). 
This indicates that the magnetic arms are much more tightly 
wound than the material arms. The effect of a finite r is 
to strengthen the non-axisymmetry, and also to slightly in- 
crease the variation of the radial phase. In addition, the 
eigenfunction extends out to somewhat larger radius in the 
t — lo/uo case than it does for r — > 0. This is illustrated by 
the second maximum with S > at r ~ 10 kpc that occurs 
when r = Iq/uq, and is caused by the tail of the spiral mag- 
netic field wrapping around an extra half-circle from where 
it crosses the corotation circle. Also shown in Fig. 9b as 
dashed lines is the envelope of the solid lines as <f> changes. 3 

The difference between the azimuthal positions where 
Btj, (or B r ) and a are maximum, denoted (or A r ), are 
shown in Fig. 9c. Positive values of these quantities imply 
that the magnetic arms lead the material arms, while neg- 
ative values imply that they trail them. The position r max 
of the global maximum of S (maximum of the dashed en- 
velope in Fig. 9b), which is the radius at which the non- 
axisymmetric component of the field is most important rel- 
ative to the axisymmetric component, is indicated with a 
vertical dotted line (red for r — >■ and black for r = Iq/uq): 
r max — 8.3 kpc for r — > and 8.7 kpc for r = lo/uo- One can 
see that r max is nearer r cor for r — > than for r = Io/uq. 

In the t -> case, studied by MS91 and SM93, it was 
found that = A r = at r = r COT and r max = r cor , 
i.e. magnetic and material arms intersect at the corotation 
radius, and magnetic arms are strongest there. Therefore, 
magnetic arms go from leading the a-spiral in the sense of 



3 The envelope is not perfectly symmetrical about 5 = because 
of interference between the m = 2 mode and the much smaller 
m = 4 mode. 
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Figure 9. Properties of the field in the kinematic regime for 
Model E at a time 3t^ g after the simulation is begun. Thick black 
illustrates the r = Iq/uq case while thin red illustrates r — > 0. In 
all plots, the corotation radius is shown as a dashed vertical line, 
a) The normalized azimuthal (solid) and radial (short dashed) 
components of the field at the azimuth <j> = (f> CO r for which the a- 
spiral crosses corotation. b) The ratio 5 of the non-axisymmetric 
to axisymmetric ^-component of the field at </> cor (solid), and at 
azimuthal cxtrcma (dashed) . The quantity (r, </> C or ) is shown 
as a thin dash-dotted line for reference, c) The phase of the (jy 
component (solid) or r-component (dashed) of the magnetic spi- 
ral minus the phase of the ct-spiral (smoothed over 5 radial grid 
points, or 100 pc). The radius r max at which the global maxi- 
mum of 5 occurs is plotted as a vertical dotted line for each r 
case, d) Negative of the pitch angle at <f> = <j> cor (solid), with 
azimuthal mean (dashed). 



Figure 10. As Fig. 9 (Model E), but now for t = lOGyr (see 
Fig. 8), when the field has reached equilibrium as a result of 
dynamical quenching. Note that in (a) The scale is now lin- 
ear instead of logarithmic. The long dashed curve is equal to 
O.lBe, 

corotation in the saturated state. 



shift. As mentioned earlier, non-locality in time produces 
a delay in the a effect of order r, so that the a-spiral has 
rotated through an angle « Qt before the dynamo has had 
the chance to respond to it. We recall that a phase shift 
of ±90° would correspond to peaks of the magnetic arms 
located mid-way between the spiral arms of a. 
The magnetic pitch angle 



the rotation for r < r cor , to trailing it for r > r COT . In the 
vanishing r case in our model, r max is slightly larger than 
r cor , and and A r are not precisely zero at r max nor at 
r cor . These differences between the results of our initial value 
treatment and those of the eigenvalue treatment of those au- 
thors can be attributed to interference from the dominant 
axisymmetric eigenmode, which is concentrated at smaller 
radius near where the dynamo number peaks, in our model. 
In any case, for both r — s> and r = Io/uq, the magnetic 
and a-spiral arms intersect in the vicinity of corotation, and 
the magnetic arms are also strongest in the vicinity of coro- 
tation, as in the earlier models. 

What is more noteworthy, however, is the phase shift in 
the magnetic arms produced by a finite dynamo relaxation 
time. This is best evaluated by comparing the values of 
and A r at the radius r max (which itself depends on r), where 
the magnetic arms are strongest. (Since B$ dominates over 
B r , A^ is more important than A r , but in the vicinity of 
the corotation radius, where non-axisymmetry is important, 
A^ ~ A r .) We find A</,(r max ) ~ A r (r max ) = —3° for r — > 
and —30° for r = Iq/uo, for an overall phase shift of —27°. 
This is comparable in magnitude to the value Qt — 18° that 
results from an order of magnitude estimate of the phase 



B r 

Pb = arctan 

Bcf, 

is plotted (by magnitude) in Fig. 9d, with its azimuthal 
mean shown as a dashed line. The negative value of ps in- 
dicates that magnetic lines have the form of a trailing spi- 
ral. The magnitude of ps is larger along the a-spiral, which 
is expected since B r is produced more efficiently where a 
is larger. The azimuthally averaged value of \pb\ decreases 
with radius such that tanps oc 1/h, as expected (Ruzmaikin 
et al. 1988; Shukurov 2005). 

In the kinematic stage, various eigenmodes are com- 
peting until the eigenfunction eventually gets smoothed out 
by the radial diffusion and grows uniformly with the global 
growth rate F. As a consequence, the excess seen near the 
corotation radius in Fig. 9a becomes less and less apparent 
with time (but returns in the saturated state, see below). 
Furthermore, the values of r max gradually increase before 
settling down (by ~ 7id,o after the simulation has begun) at 
nnax = 9.3 kpc for r — > and r max = 9.5 kpc for r = Iq/uq. 
The key results mentioned, in particular the large phase shift 
caused by a finite dynamo relaxation time, are valid for any 
time during the kinematic regime. 
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5.2 Magnetic arms in the nonlinear dynamo 

Many of the features of the kinematic solution are also ap- 
parent in the nonlinear regime. Figure 10 is similar to Fig. 9, 
but now for the saturated dynamo. In Fig. 10a, it can be seen 
that nonlinear effects significantly enhance the m = 2 part 
of the mean magnetic field near the corotation radius. In 
Fig. 10b, it is clear that deviations from axisymmetry are 
only important within « ±2kpc of the corotation radius, 
and that the corresponding eigenfunction is strongly peaked 
about the corotation circle. Fig. 10c shows the phase shifts 
and A r , with the location of the global maximum r max 
of 8 shown as a dashed vertical line (red for r — > and 
black for finite r). In this case, r max ~ 8.0 kpc for r — > 
(consistent with the corotation radius) and r max ~ 8.6 kpc 
for t — lo/uo- Because the eigenfunction is sharply peaked, 
it is more useful to consider the values of and A r at 
r max rather than their asymptotic values far from corotation. 
The phase differences A r (r max ) and A^,(r max ) are +3° when 
r — > and —27° for r = lo/uo Thus, the non-axisymmetric 
parts of the magnetic field are substantially displaced from 
the Q-spiral arms by an overall phase shift of —30° due to 
the finite value of r. The magnitude of the displacement is 
again comparable to fir = 18°, as expected. Fig. lOd shows 
the magnitude of the pitch angle at 4> — ^or as well as the 
azimuthal mean of this quantity for both r cases. When com- 
paring to the dash-dotted curve illustrating a k at (j> = 4> COT 
in panel (b), it is clear that, like in the kinematic regime, 
the magnitude of the pitch angle is generally larger along 
the Q-spiral. However, the decreasing trend ps oc l//i is no 
longer present in the nonlinear regime, and the azimuthal 
average of ps (dashed curve) is almost constant over a large 
range of radius for both r cases. 

Figure 11 shows the magnetic field strength (upper pan- 
els), the quantity 8 (middle panels), and pitch angle (lower 
panels) obtained in Model E, with the case of r — > on 
the left, r = Iq/uq in the middle column, and r = 2?o/«o 
on the right. The upper panels contain information similar 
to that of Fig. 10a, the middle row of panels are similar to 
Fig. 10b, while the lower panels are similar to Fig. lOd. In 
the upper panels, which show the field magnitude, magnetic 
arms are clearly visible near the corotation circle. As they 
are more tightly wound than the a-arms, they cut across 
them, as noted above. Comparing panels (a)-(c) or (d)-(f) 
of Fig. 11, and also the black and red curves of Fig. 10b, 
we may note the changes to the magnetic field caused by 
a finite r, particularly outside of the corotation circle. The 
magnetic arms are more pronounced (larger B at their cen- 
tres) , are more sharply defined (faster variation of the radial 
phase), and can be seen to extend for a longer azimuthal 
angle outside of the corotation circle for the finite r case 
compared to the r — > case. Since the mean electromotive 
force responds to an enhancement in a with a time delay r, 
the maxima in magnetic field strength are displaced down- 
stream from those of a by an angle of order fir. Moreover, in 
the r = lo /uq case, a more significant part of the magnetic 
arm clearly lies in between the a-spiral arms. These features 
can be seen more clearly for r = 2/o/wo in the right hand 
panels of Fig. 11. As r is increased, the downstream (trailing 
the a-spiral) part of the magnetic arms becomes more and 
more extended, whereas the upstream part is suppressed. 

The magnitude of the pitch angle is displayed in the 



lower panels of the figures, and clearly also has a spiral pat- 
tern. The pitch angle, as opposed to the field magnitude, 
does not vary strongly with radius. Because of this, the spi- 
ral morphology of the field is more clearly visible in the 
pitch angle. However, an accuracy of a few degrees would be 
required to measure the differences. 

Another quantity that can be obtained is the pitch angle 
of the magnetic ridges, i.e. of the magnetic arms themselves 
as opposed to the magnetic field that constitutes them. 
This can then be compared with the pitch angle of the a- 
spiral arms. To estimate the pitch angle of the a-spiral, from 
Eq. (27), we may write down an expression for the values 
of azimuth (along one of the n arms) <^ max which maximize 
ak, 



hftt. 

n 



(39) 



The pitch angle p a of the alpha spiral can then be obtained 
using the expression 



cot p 



dr 



(40) 



Substituting the values from Model E, k = —20R~, n = 2, 
as well as r = r cor = 0.4R, we obtain p a {r CO i) = —14°. 

Now we obtain the pitch angle of the magnetic ridges. 
We define a ridge as the spiral connecting the maxima of the 
magnetic energy E, or, equivalently, the field magnitude. 
Using the same approach as in Eq. (40) above, the pitch 
angle ps of one such ridge can be written as 

dA E 



cot pe = — + r- „ 
n or 



(41) 



where the phase shift Ae of the ridge with respect to the 
a-spiral is analogous to the quantities A r and A^ discussed 
above. We evaluate the pitch angles at r = r max , where 
the global maximum of the quantity 8 is located, recalling 
that r max = 8.0 kpc for r — > 0, while r max = 8.6 kpc for 
r = lo/uo- We then find jJ.e(r max ) ~ —5° for r — ► 0, and 
PB(r max ) ~ — 4° for t = Zo/mo- The pitch angles of the mag- 
netic ridges are quite small, even as compared to the rather 
tightly wound a-spiral. One must keep in mind, however, 
that the values obtained for pe are sensitive to the shear 
near corotation. 



5.3 Coupling of the spiral pattern to galactic 
dynamos: multitude of mechanisms 

Modulation of the a effect is just one of the mechanisms 
through which the spiral pattern can affect the dynamo ac- 
tion. In this section we consider two other plausible mecha- 
nisms in more detail, which eventually have to be combined 
into a conclusive picture. 

5.3.1 Outflow speed modulated by the spiral pattern 

Model F has the vertical advection velocity enhanced along 
a spiral (with the a-spiral turned off), to model stronger 
galactic outflows from the spiral arm regions, where star for- 
mation is enhanced. We use Ru — 0.90 at r = R/2 = 10 kpc 
and eu = 0.5 so that the outflow speed within the arms, 
corresponding to Ru = 1.35, is too large for the optimal 
field growth near corotation, while the inter-arm value of 
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Figure 11. The top panels (a)-(c) show the magnitude of the magnetic field in the saturated state (t = lOGyr) in the r — > (left), 
t = Iq/uq (middle) and r = 21q/uq (right) cases. In all panels, the peak of c% is shown as a thick black line, while the corotation radius 
is illustrated as a black dotted circle. The colour of the central regions in (a)— (c) has been saturated in order to visually bring out the 
behaviour near corotation, where the non-axisymmetric modes are important. Panels (d)-(f) show the ratio 5 of the non-axisymmctric 
to axisymmetric components of B$. It can be seen that the r effect acts to enhance the non-axisymmetry in the field. It also causes a 
phase shift of the magnetic arms toward the trailing side of the material arms, leading to magnetic arms that are more pronounced in 
the region between the material arms than for the r — > case. The bottom panels show the magnetic pitch angle. Again, the colour 
within a region around the center has been saturated to render the non-axisymmetry of the pitch angle more apparent. 




B/B 



0.00 0.02 0.04 0.06 0.08 0.10 

Figure 12. As Fig. 5, but now for Model G (random seed with spiral forcing). The r — ¥ case is shown. The peak of is shown as a 
thick black line, while the corotation radius is marked by the black dotted circle. In addition to reversals, strong deviations from axial 
symmetry develop near the corotation radius. The morphology of the magnetic field is affected by the outward-moving reversals. 
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Figure 13. Evolution of the magnitude of the magnetic field in Model H (zero helicity flux, U z = 0) and r — > 0: (a) t = 2.2 Gyr, 
(b) 3.3 Gyr, (c) 4.4 Gyr, and (d) 5.5 Gyr. The peak of 0% is shown as a thick white line, while the corotation circle is indicated with 
a white dotted circle. The expanding ring of the mean magnetic field acquires a large m = 2 component as it passes by the corotation 
circle. 



Ru = 0.45 is close to being optimal. This mechanism for 
producing magnetic arms located between the material spi- 
ral arms was first suggested by Sur et al. (2007). 

We find that this mechanism can indeed lead to a non- 
axisymmetric field which peaks in the inter-arm regions. For 
the region of parameter space explored, the effect is, how- 
ever, too weak to be of much consequence, though it might 
still be important for different parameter values. 

5.3.2 Turbulence modulated by the spiral pattern 

The characteristic strength of magnetic field at which non- 
linear dynamo effects become pronounced, B eq , can also be 
affected by the spiral pattern, e.g., through the variation in 
the turbulent energy density. We considered a model with 
B eq enhanced along a spiral in the same way as above for 
and U z . Since B cq = (Anp^^u, and we keep u constant, this 
is tantamount to modulating p. Although an enslaved spiral 
does result, we find the effect to be too weak to be of much 
consequence, at least for the parameter space explored. 

Other possibilities have been explored in the literature, 
e.g. the modulation of the turbulent magnetic diffusivity %, 
but a more extensive study of these effects is left for a future 
work. 

5.4 Exploring the parameter space 

The purpose of Models G-L is to test the results under 
reasonable variations of parameters. We have varied the 
strength of the seed field, the mean vertical velocity U z , 
the strength and the pitch angle of the spiral pattern via e a 
and k, respectively, and the number of material spiral arms 
n. 

5.4.1 The seed field 

In Model G, we use a seed field that is random Gaussian 
noise. Model G differs from Model C in that the underlying 
disc in non-axisymmetric. We find reversals to occur in more 
or less the same locations as in that model. However, as can 
be seen in Fig. 12 (r — > case), a strong quadrisymmetric 
magnetic field component is now apparent, superimposed 



on an axisymmetric pattern with reversals. The morphology 
of the magnetic spiral arms is significantly affected by the 
reversals, as the horizontal magnetic field is necessarily zero 
wherever a reversal occurs (compare, e.g., with Fig. 11a for 
Model E). 

5.4-2 Gas outflow and magnetic helicity flux, U z 

The effect of vanishing helicity flux is discussed in Sect. 4.1.1 
for the axisymmetric dynamo (Model A). Model H is de- 
signed to consider this for non-axisymmetric magnetic fields, 
and thus differs from Model E only in that it has U z = 0. 
This situation can occur in galaxies after the end of a burst 
of star formation, when the galactic fountain or wind ceases 
and the mean-field dynamo action is choked by the mag- 
netic helicity conservation to leave magnetic field decaying, 
assuming no other flux is important. 

As expected, the mean field decays after a period of tem- 
porary growth, as shown in Fig. 13 for the nonlinear regime. 
The expanding annular magnetic structure of Sect. 4.1.1 is 
prominent, but develops strong deviations from axial sym- 
metry as the ring approaches and then passes through the 
corotation circle. At late times (not shown), the magnetic 
ring fades in intensity by about a factor of two in the 
next 6 Gyr) and becomes more and more axisymmetric as it 
moves out. 

5.4-3 The strength of the spiral pattern, e a 

In Model I, e a is changed from 0.5, as in Model E, to e a = 1, 
so that now — between the spiral arms. The saturated 
solution is illustrated in Fig. 14. As could be expected, this 
enhances the deviation of the mean magnetic field from axial 
symmetry, manifested in more pronounced magnetic arms 
and a larger azimuthal variation in the magnetic pitch angle. 

The magnitude of the non-axisymmetric part of the 
magnetic field mode slightly exceeds the axisymmetric part 
near r = r max when r = Iq/uq. This produces small regions 
of positive pitch angle just outside of the corotation, visible 
in Fig. 14f (darkest regions), where the magnetic lines locally 
have the shape of a leading spiral. Putting r = Iq/uq also 
results in a somewhat larger phase shift compared to the 
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Figure 14. Same as in Fig. 11, but now for Model I (e a = 1 
instead of 0.5), and with the rightmost column of panels omitted. 
The first column of panels (a, c, e) is for r — ¥ while the second 
column (b, d, f) is for r = lo/uo- The non-axisymmetry of the 
field is stronger than it is for t a = 0.5. Note that the very darkest 
parts of the blue regions in (f) have positive magnetic pitch angle, 
i.e., magnetic lines in the form of a leading spiral locally. 

r — lo/uo case in Model E. At the radius where the devia- 
tions from axial symmetry are the strongest, r = r max , we 
obtain A r (r max ) ~ A^(r max ) = -36° for r = l /u Q , whereas 
for Model E we had -27°. On the other hand, for r -> 0, 
the phase differences are both equal to +3° , the same as in 
Model E so the overall phase shift is —39°. The trailing part 
of the magnetic spiral is clearly enhanced when r is finite. 

5.4-4 The pitch angle of the material spirals, k 

In Models J and K we make the a-spiral less tightly wound 
by reducing the parameter k. It is evident from Fig. 15 that 
changing the value of k from — 20i? _1 (p Q (r cor ) = —14°, 
Model E) to -8i? _1 (p a (r cor ) = -32°, Model J) does not 
have any significant qualitative effect on the magnetic field. 
Even replacing the material spiral with a bar (ft = 0, 
Model K) does not lead to a drastic change in the satu- 
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Figure 15. Same as in Fig. 14 but now for Model J (e a = 0.5 
and k = — 8.R instead of —20R~ 1 ). Magnetic arms are strong 
in the regions in between the material arms (t -> on the left 
and T = lo/uo on the right) . 

rated magnetic field, so that bars also lead to spiral mag- 
netic fields. This happens because of the differential rotation 
of the gas, which shears out the enhancement of the field due 
to the bar. 

The relative strength of the non-axisymmetric part of 
the magnetic field is slightly larger in a disc with a more 
open spiral pattern, which could be expected since the ra- 
dial magnetic diffusion becomes weaker. An important effect 
of a more open spiral pattern, visible in Fig. 15a-b, is that 
magnetic spiral arms are now mostly in between the mate- 
rial arms, although still confined to an annular region around 
corotation. The values of A r (r ma x) and A^(r ma x) are both 
—33° for t = lo/uo, and +3° for r — > 0. Thus, the finite re- 
laxation time of the mean emf leads to a significantly larger 
phase shift than in Model E (—36° vs. —27°). The variation 
of the magnetic pitch angle is strongly modified, because 
the pitch angle is enhanced where a is large. However, the 
range of variation in the magnetic pitch angle remains al- 
most unchanged from the case of the more tightly wound 
spiral pattern of Model E. 
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Figure 16. Same as in Fig. 14 but now for Model L (e a = 0.5 
and n = 4 instead of 2). Magnetic arms are again strong in the 
regions in between the material arms. For the t = Iq/uq case 
(right panels), the magnetic arms are more pronounced than for 
r — > 0, and are mostly trailing the corresponding material arms. 

5.4-5 The number of spiral arms, n 

We tried various values of n, both odd and even, and found, 
unsurprisingly, that the number of magnetic arms is equal to 
the number of material arms. The results of Model L (n = 4) 
are plotted in Fig. 16. The magnetic arms are largely located 
in between the material arms, especially for r = Iq/uq. They 
are also much stronger and more well-defined in the finite r 
case. 



6 MAGNETIC RESPONSE TO TRANSIENT 
SPIRAL PATTERNS 

6.1 A rigidly rotating spiral pattern 

As galactic spiral patterns may be transient in nature, we 
now explore, in Model M, how the magnetic field responds 
to sudden changes in the spiral forcing. The runs are identi- 
cal to Model B (axisymmetric disc) up until the a-spiral is 
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Figure 17. As Fig. 8, but for Model M (transient rigidly rotating 
spiral). The times at which the a-spiral was turned on and off are 
indicated by vertical dashed lines. 



turned on, at the time t on = 5.1 Gyr, when the dynamo is 
already in its nonlinear phase. Subsequently, the parameters 
of the model are identical to those of Model E (standard spi- 
ral forcing) up until the time t a s = 7.3 Gyr, when the spiral 
modulation of a is turned off. 

The evolution of the magnetic field strength in each 
Fourier mode, averaged over the area of the disc, is shown 
in Fig. 17. The m = 2 mode, which initially is present solely 
because of numerical noise, responds rapidly to the onset of 
the a-spiral, and then behaves almost exactly as in Model E 
(compare with Fig. 8). In fact, for r = lo/uo, magnetic en- 
ergy in the m — 2 part of the mean field slightly exceeds 
that in Model E before being reduced to the latter. The 
timescale of the adjustment of the magnetic field, i.e., the 
time from the onset of the a-spiral, taken for the m — 2 part 
in Model M to grow up to that in Model E, is as short as 
~ 0.2 Gyr. This timescale is shorter than the vertical turbu- 
lent diffusion time, but comparable to the dynamo growth 
time at corotation. 

Figure 18 shows characteristic magnetic configurations 
for t = lo/uo (the r — >■ case is qualitatively similar so 
only r = lo/uo is shown). After a short period of 0.15 Gyr 
after the onset of the spiral forcing, the deviation from ax- 
ial symmetry in the magnetic field is strongest somewhat 
inside of the corotation circle, as can be seen by the lo- 
cation of the peak in amplitude of the dashed envelope in 
Fig. 18a. This is explained by the dynamo responding more 
rapidly to sudden changes in the disc at locations where 
the dynamo number is larger. At t = t on + 0.29 Gyr, the 
non-axisymmetric component of the field has expanded out- 
ward in radius, and has become stronger as compared to 
the axisymmetric component (panels b and f). This is to 
be expected since, according to results already discussed, 
deviations from axisymmetry are most important near the 
corotation radius. Next, at t = ton + 2.2 Gyr = t s (Fig. 18c 
and g), the magnetic field has reached nearly the same equi- 
librium state as in Model E (compare with Figs. 10b and 
lib). After the a-spiral has been switched-off, the magnetic 
arms remain for a rather long time compared with the time 
it had taken for them to arise. This can be seen in Figs. 18d 
and h, where we show the field at t = t s + 0.37 Gyr. The 
energy in the m — 2 part declines to a quarter of its value at 
toff within 0.2 Gyr and 0.3 Gyr, respectively, for r — > and 
r = Iq/uq. Therefore, 'ghost' magnetic spiral arms remain 
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Figure 18. Evolution of the magnetic field under the action of a transient, rigidly rotating o-spiral (Model M) for r = Io/uq. Each 
column shows the ratio of non-axisymmetric to axisymmetric part of B$, as in Fig. 9 (top) and magnitude of the field as in Fig. 11 
(bottom) at the following times after the emergence of the ct-spiral: (a, e) 0.15 Gyr, (b, f) 0.29 Gyr, (c, g) 2.2 Gyr, followed by (d, 
h) 0.37 Gyr after the a-spiral is turned off. 



after the demise of the material spiral arms, and survive 
longer when r is finite. The extra delay of about 0.1 Gyr 
due to the finite r is larger than the delay r ~ 0.01 Gyr that 
might naively be expected. 

In summary, (i) a magnetic spiral can arise within a 
short period of 0.1-0.2 Gyr after the sudden emergence of 
an a-spiral, (ii) the magnetic spiral may still be obvious a 
few hundred Myr after its turnoff, especially if r is large, and 
(iii) magnetic field responds to the onset of a transient ma- 
terial spiral first inside of the corotation circle, then further 
outward near the corotation, and then beyond. 



6.2 A winding-up spiral pattern 

The models described so far had a rigidly rotating spiral 
with a constant pattern speed. As a consequence, the mag- 
netic spiral arms are localized within a few kpc of the coro- 
tation radius. However, magnetic arms, e.g., in the galaxy 
M51, can be coherent from within about 1 kpc of the galac- 
tic centre, out to 10 kpc (Fletcher et al. 2011). Observations 
of M51 indicate that there does not exist a single pattern 
speed for the spiral (Shetty et al. 2007), and simulations 
of that galaxy suggest that the spiral arms wind up with 
time (Dobbs et al. 2010). Some N-body simulations of disc 
galaxies find that material spirals have 'pattern speeds' de- 
creasing with radius. In some cases, the radial profile of the 
rotational velocity of the spiral almost matches the rotation 
curve of the galaxy (Grand et al. 2011). The rigidly rotat- 
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Figure 19. Similar to Fig. 10 (with changes to the plotting 
ranges), but now for Model N (the ojj spiral winding up with 
the gas). The figure shows the field at 0.1td,o — 73 Myr after the 
sudden onset of an 0% bar when the dynamo action has already 
saturated. 
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Figure 20. Same as in Fig. 19, but at 0.2t d ,o =! 146 Myr after the 
onset of the bar, (which subsequently winds up into a spiral). 
The non-axisymmetry is enhanced from the earlier time, and now 
extends over a larger radial range. 



ing spiral model may best characterize some galaxies, while 
models with spiral winding up with the differential rotation 
may be more suitable for others - in any case they plausi- 
bly represent two extreme cases. A winding-up spiral may 
also be thought of as the large- N limit of N non-overlapping 
spiral patterns in the disc. 

We have already explored rigidly rotating a-spirals, so 
we now turn to the case of an a-spiral that winds up with 
the gas. In this model every radius effectively becomes a 
corotation radius, so based on our results for a rigidly ro- 
tating spiral, we would expect the non-axisymmetric part 
of the magnetic field to be more extended radially, as seen 
in some observations. In Model N, the field evolves as in 
Model B (axisymmetric disc) until the time f on , when we 
use parameter values e a = 0.5 and n = 2, as in Model E, 
but now let the spiral wind up starting from the 'bar' stage 
(ft = 0) at t = t on = 6.6 Gyr, when the m = field is already 
saturated. We find that strong deviations from axial sym- 
metry in the magnetic field develop during the first 150 Myr 
after t on . In Figs. 19 and 20 we show the properties of the 
resulting field at two times, t = t on + O.ltd.o = t on + 73 Myr 
and ton + 0.2£d,o = ton + 146 Myr, respectively. 

As can be seen by comparing panels a and b of Fig. 19 
(t = t n + 73 Myr) and also by comparing panels a and b 
of Fig. 20 (t — ton + 146 Myr), the non-axisymmetric and 
axisymmetric parts of the magnetic field are the strongest 
at about the same radii, so that the magnetic field is es- 
sentially non-axisymmetric. This can also be seen from the 
colour plots of Fig. 21 (t = t on + 73 Myr) and Fig. 22 
(t — t n + 146 Myr). In addition, the magnetic spiral arms 
are extended throughout the disc over a much larger range 
of radii than with a rigidly rotating spiral pattern (where the 
non-axisymmetric field concentrates around the corotation 
circle). 



Figures 19b and 20b show that the radial phase of the 
field (solid lines) is similar to that of ak (dash-dotted line), 
and that the two are well-correlated, demonstrating that the 
dynamo responds rather quickly to the spiral forcing. From 
Figs. 19c and 20c, we see that for r — > the phase difference 
between the magnetic and material spiral arms is negligible 
at radii where the magnetic spiral arms are strong (near 
where the red solid and dashed lines cross the red dotted 
vertical line). However, for finite r = lo/uo, the magnetic 
arms lag the material ones by a large angle, A r ~ A$ ~ 
— 15° to —25° over a range of radius of order 10 kpc. The 
values vary somewhat with radius and with time. This phase 
shift is strikingly apparent when comparing panels (a) and 

(e) of Fig 21 for r — s> with, respectively, panels (b) and 

(f) for t — lo/un. Only in panels (b) and (f) do the yellow- 
orange regions of high field strength lag the black centers of 
the cv-spiral arms by a significant amount. The same is true 
in Fig. 22, showing snapshots at t = t on + 146 Myr. 

As might be expected, the phase shift is larger for r = 
21q/uo (—20° to —30°) and even larger for r = Mo/uo (—25° 
to —35°), as can be seen from the third and fourth columns 
of Figs. 21 and 22. However, as r is increased, the non- 
axisymmetric part of the magnetic field tends to weaken, as 
the spreading of the response of the dynamo in time due 
to the finite dynamo relaxation time is translated into an 
effective blurring of the a-spiral in space (the effect is most 
evident when comparing panels e-h of Fig. 21 or Fig. 22). 

We also tried, separately, e Q = 1 and n — 4. For e a = 
1, the results are similar to Model N (e a = 0.5), with a 
comparable magnitude of the phase shift in the r = Iq/uq 
case. However, we find magnetic arms that are somewhat 
stronger and more extended than in Model N, as might be 
expected from the discussion of Sect. 5.4.3. On the other 
hand, when n is increased from two to four, we find a four- 
arm magnetic spiral, as might be expected. Importantly, due 
to the smaller angle between the material arms, the magnetic 
arms are located closer to the centres of the inter-arm regions 
than for the two-arm spiral of Model N (see also Sect. 5.4.5) 

This model, which incorporates a spiral pattern that 
seems to wind up with the differentially rotating gas, along 
with a dynamo relaxation time r > lo/uo, seems to be the 
most successful, among the models explored in the present 
work, at explaining the large phase shifts between mag- 
netic and material arms, roughly constant over a large radial 
range, observed in some galaxies, such as NGC 6946. 



7 DISCUSSION AND CONCLUSIONS 

We have presented a set of models of the mean-field galac- 
tic dynamo action that include: (i) dynamo forcing by sta- 
tionary and transient spiral structure, (ii) time delay in the 
response of the mean electromotive force to changes in the 
mean magnetic field (finite relaxation time r), and (iii) dy- 
namo nonlinearity based on magnetic helicity balance. 

We incorporate into a detailed galactic model a phys- 
ically and mathematically justifiable dynamo saturation 
mechanism based on the advection of magnetic helicity 
by galactic outflows (Shukurov et al. 2006; Sur et al. 
2007). Within this framework, we explore the behaviour 
of the mean magnetic field in both axisymmetric and non- 
axisymmetric galactic discs. A novel aspect of this paper is 
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Figure 21. Winding-up spiral of Model N. The figure is similar to Fig. 11. From left to right the columns represent: r — > 0, r = Io/uq, 
t = 21q/uo, and r = 4Zo/«o- All panels are snapshots at O.ltd.o = 73Myr after the sudden onset of the 'bar'. Increasingly larger 
phase shifts in the magnetic arms toward the inside of the spiral arms (trailing side) are evident as r increases. 



the introduction of a finite relaxation time into the galactic 
mean-field theory and a detailed exploration of its effects, 
particularly on the non-axisymmetric dynamo modes. We 
discuss how the imprints of the galactic spiral pattern on 
the morphology of the mean magnetic field are affected by 
the finite relaxation time of the dynamo. Unlike many ear- 
lier studies of galactic dynamos, we approach these problems 
from the evolutionary viewpoint and consider the magnetic 
response to both steady and transient spiral patterns, and 
briefly discuss the diverse ways in which galactic spiral arms 
can affect the dynamo action. 

We have confirmed, in the new framework, earlier re- 
sults on the evolution of axisymmetric magnetic fields in ax- 
isymmetric discs, including the spreading of magnetic fronts 
along radius and the occurrence of reversals. As expected, 
the magnetic field is able to spread radially to encompass the 
whole disc within the galactic lifetime, as long as the mag- 
netic helicity flux is sufficient, but not too large. What is new 
is that in the absence of any helicity flux, contracting and 
expanding rings of magnetic field, emanating from the re- 
gion where the magnitude of the dynamo number is largest, 
can propagate for several Gyr in the nonlinear regime. This 
is due to the variation along radius of the time taken for the 
magnetic field to reach saturation. 



The lifetimes of reversals in the galactic disc have been 
known to be sensitive to the details of the galactic rota- 
tion curve, the geometric shape of the gas layer, details of 
the turbulent energy distribution in the galaxy, etc. We have 
added another important parameter to this list, the dynamo 
relaxation time. We also find that such reversals can have 
non-trivial effects on the morphology of non-axisymmetric 
magnetic structures in a non-axisymmetric disc. It would 
be interesting to carry out a more detailed study of mag- 
netic fronts and the possibility of long-lived reversals in the 
context of the above framework. 

We have also shown that it is difficult to maintain non- 
axisymmetric magnetic structures in an axisymmetric disc, 
which leaves one little choice but to try to explain such struc- 
tures as arising from a non-axisymmetric disc. Most of our 
effort was then directed to the response of the mean mag- 
netic field to the galactic spiral pattern. The focus in this 
part of our work is the nature of the so-called magnetic 
arms, spiral-shaped regions of enhanced regular magnetic 
field (traced by polarised radio emission and, perhaps, Fara- 
day rotation) that are located in between the gaseous and 
stellar spiral arms in some galaxies (e.g., NGC 6946 and 
IC 342) but can overlap with or cross the material arms in 
other galaxies (e.g., M51). We model the effect of the spiral 
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Figure 22. As Fig. 21 but now for t = t on + 0.24^ = *on + 146 Myr. In the middle and bottom rows, the colour in the central region 
has been saturated to bring out the behaviour at larger radius. 

pattern on the dynamo by assuming that the a effect or the 
galactic/fountain outflow or the equipartition field are en- 
hanced within the spiral arms. A non-axisymmetric a effect, 
in particular, does lead to strong magnetic arms. 

When the dynamo relaxation time r — > 0, the mag- 
netic and material arms almost coincide at corotation, as 
expected from earlier work (MS91; SM93; M96). However, a 
finite relaxation time causes a significant azimuthal lag be- 
tween the magnetic and material arms at the radius where 
the magnetic arms are strongest (near corotation). In this 
respect, our model can explain a wider range of observa- 
tions than earlier models. However, we concur with the ear- 
lier authors in that the non-axisymmetric parts of the mean 
magnetic field driven by the spiral pattern are mainly lo- 
calised around the corotation radius, extending about 4 kpc 
in radius. The radial extent is, of course, model-dependent 
and can be larger in regions with weaker rotational velocity 
shear. 

The corotation radius is also approximately where each 
magnetic arm crosses the corresponding material arm, going 
from leading (with respect to the direction of the galactic 
rotation) the material arm inside the corotation circle to 
trailing it outside the corotation. The primary effect of the 
finite relaxation time, in the case of a steady, rigidly rotating 
spiral pattern, is to suppress the leading part of the magnetic 



spiral arms and to enhance, and extend azimuthally, the 
trailing part. Thus, this effect produces a prominent 'tail' 
of large-scale magnetic field in between the material arms. 
This is even more true when the number of arms is increased 
from two, to say, four. 

The magnetic arms of the spiral galaxy NGC 6946 are 
'interlaced with' (in between) the material arms (Beck & Ho- 
ernes 1996; Beck 2007), and have been called 'phase-shifted 
images' of the preceding (in the sense of the rotation) ma- 
terial arms (Frick et al. 2000). Although the implied phase 
shifts are more constant with radius than in our model, we 
have obtained a shift of approximately the same magnitude 
(30°-40°) and in the right direction. More generally, the 
shift is of order Qt, where Q is the pattern angular velocity 
of the material spiral. Larger values of r thus lead to larger 
phase shifts. This is caused by the delay in the a effect of 
order r, so that by the time the dynamo has had the chance 
to respond to the enhancement in a along the material arm, 
the arm has already rotated by an angle ~ fir. 

Spiral density waves may be transient and we con- 
sider how quickly the regular magnetic field can respond 
to changes in the galactic spiral structure. We find that the 
response time is small (< 0.2 Gyr). On the other hand, we 
find that magnetic spiral arms survive for several hundred 
Myr following the destruction of the material spiral arms 
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in the dynamo nonlinear regime. Moreover, a finite dynamo 
relaxation time is found to significantly prolong the life of 
such lingering magnetic arms. This opens the intriguing pos- 
sibility of 'ghost' magnetic arms, which were produced by 
material arms that have since disappeared. 

Despite the wide range of models considered, our suc- 
cess in reproducing magnetic arms interlaced with the mate- 
rial arms as perfectly as it is believed to happen in NGC 6946 
is admittedly limited in models assuming forcing of the dy- 
namo by a steady, rigidly rotating spiral. In addition, this 
type of spiral forcing leads to magnetic arms concentrated 
over a smaller range in radius than is observed in many 
galaxies. Another possibility, rendered more likely by sev- 
eral recent studies (e.g. Dobbs et al. 2010; Sellwood 2011; 
Quillen et al. 2011), is that the spiral patterns of many 
galaxies, rather than being rigidly rotating, as is usually as- 
sumed in galactic dynamo models, in fact wind up (at least 
to some extent). This may happen if, for instance, there 
are interfering two and three-armed spirals rotating at dif- 
ferent angular frequencies (A. Quillen, private communica- 
tion) . The two-arm grand design spiral pattern of the galaxy 
M51 (NGC 5194), thought to be caused by tidal forcing by 
its neighbour, NGC 5195, is also found to wind up in detailed 
N-body simulations that are able to accurately reproduce its 
spiral morphology (Dobbs et al. 2010). 

With these recent advances in spiral structure the- 
ory in mind, we also investigated the opposite extreme to 
rigidly rotating patterns: material arms that are wound up 
by the galactic differential rotation. The nonlinear dynamo 
responds very quickly to such forcing, and for vanishing 
dynamo relaxation time the mean magnetic field more or 
less traces the spiral arms over a large range in radius (as 
seen in many observations) and winds up with them. On 
the other hand, for finite dynamo relaxation time there is a 
large azimuthal lag of each magnetic spiral arm compared 
to the corresponding material arm over a large range in ra- 
dius. Magnetic arms trail the material arms by 15°-25° (for 
r = io/tio), varying somewhat with time and radius over the 
disc. This shift is of order tor, where uj is the angular velocity 
of the gas, and we have shown that larger, but still plausi- 
ble, values of r, lead to even larger phase shifts. Increasing 
the number of spiral arms also causes the (equal number 
of) magnetic arms to be located closer to the centres of the 
inter-arm regions, so that they may be described as inter- 
laced with the material arms. Therefore, allowing for the 
possibility of the spiral winding up can drastically improve 
agreement with observations of the regular magnetic fields 
in some galaxies, but with the trade-off that the magnetic 
arms (that in this model, either trace or interlace the ma- 
terial arms) are almost as short-lived as the spiral patterns 
that presumably generate them. 

In summary, we have incorporated into galactic dynamo 
theory several well-studied physical effects not previously 
considered, namely (i) non-locality in time and (ii) forcing by 
both spiral arms which steadily rotate and those which wind 
up due to differential rotation. This allows several observed 
features of magnetic arms to be more naturally reproduced. 
Particularly interesting are the models we have presented 
in which magnetic arms extend over a large range of radii 
and either trace material arms over several kpc, or else are 
phase-shifted images of material arms, trailing them in the 
sense of the galactic rotation. It would be interesting and 



important to test these ideas by applying them to specific 
galaxies, for which the rotation curve, velocity dispersion, 
spiral structure, etc., can be constrained by existing data 
from observation and simulation. 
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APPENDIX A: A MORE GENERAL FORM OF 
THE BASIC EQUATIONS 

Al The telegraph equation 

Having defined 

J== V x £, 

and assuming that r], ijt, t and c T are constants, the cylin- 
drical polar components of Eq. (3) can be written as 
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Taking the curl of both sides of Eq. (9), we get 
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A2 Dynamical quenching by advective flux 



Assuming that a — a k + a m , that is, the a effect can be 
separated into kinetic and magnetic contributions, we have 
(Shukurov et al. 2006; Subramanian & Brandenburg 2006) 
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where B^q = Aivpu 2 is the equipartition field strength and 
Ti m = rjt,o/rj is the magnetic Reynolds number. To evaluate 
the term containing £, we need the evolution equations for 
its components. From Eq. (9), we find 
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The evolution equations for £ are redundant given those for 
T ', since the latter could be obtained by taking the curl of 
the former. However, we find it convenient to solve sepa- 
rately for £ and T (see Sect. 3.1). 



APPENDIX B: THE NO-z APPROXIMATION 
Bl Vertical diffusion 

Under the no-z approximation, the second derivatives 
with respect to z can be approximated as d 2 Bi/dz 2 ~ 
—-iv 2 Bi/4:h 2 (where i — r,(f>), which gives the correct sign 
of the diffusion term. This approximation can be derived 
from the one-dimensional eigenfunctions obtained from the 
perturbation theory (Sur et al. 2007; Shukurov & Sokoloff 
2008). 



B2 The a effect 

When applying the no-z approximation to the terms 
d(aBi)/dz in Eqs. (A4) and (A5), one must be careful about 
the sign. The sign must be chosen so that the a effect can 
contribute to a positive dynamo growth rate. This logic leads 
to the adoption of the approximations 
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The second of these is only relevant when the a 2 effect is 
taken into consideration. 



Furthermore, Phillips (2001) has shown that the no-z 
approximation can be made more accurate with the addi- 
tional numerical factor 2/n, at least for the aui dynamo. As 
a matter of symmetry, we include the same numerical factor 
in front of both a terms, so that we may write 
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It may be asked whether extending the no-z approxi- 
mation to include the a term in the evolution equation for 
B$ (through Eq. A5), as we have done here for the first 
time, actually helps to improve the accuracy of the solu- 
tion. To answer this, we compared the solutions obtained 
using simple two-dimensional (in r — 2) acj and a 2 uj galactic 
dynamo model with those obtained from the correspond- 
ing one-dimensional (in r, with no-z) models. Interestingly, 
we found significantly better agreement between the ID and 
2D solutions when the a 2 effect was included in both models 
than when it was left out. 



B3 Vertical advection 

The terms — BV • U (divergence) and — U • VB (advection) 
are approximated as 
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For —d(a m U z )/dz in (A7), we have 
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B4 Approximation for £ ■ B 

For r — > 0, we have 

£ B = c r a(B 2 r + Bl) + c rVt (s r ^f - ^^7) - (B4) 

and the second bracketed term arising from the mean current 
helicity vanishes in the no-z approximation. Therefore, a 
more precise method must be used to estimate this term. 
With this in mind, Sur et al. (2007) substituted the one- 
dimensional perturbation solution of the dynamo equations, 
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and its derivatives with respect to z into Eq. (B4). The re- 
sulting expression for £ ■ B was then averaged over < 
z ^ ft. In this way they obtained a non-zero estimate for 
c T rjt(B r dB^,/dz — B^dBr/dz), which stems from the term 
of order J~^D in Eq. (B5). 

In the general case of finite r, the correction must come 
in the evolution equations (A8) and (A9) for £ r and 
Averaging Eqs. (B5) and (B6) and also their z-derivatives 
over ^ z ^ ft, we find 
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dB 4 
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These expressions are used in Eqs. (A8) and (A9) for calcu- 
lating £ ■ B in Eq. (A7). 

Substituting expressions (B7) and (B8) into Eq. (B4), 
we find for the r — > limit, 
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'8^/2/! 

The term involving rj t used here is smaller by a factor -k 
than that of Sur et al. (2007) since we have taken the scalar 
product of the 2-averages of £ and B, whereas those authors 
used the average of the scalar product. This difference does 
not appear to be important. 
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